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I . INTRODUCTION 

The  Linus  program  at  NRL  has  as  its  goal  the  design  and 

demonstration  of  a fusion  reactor  utilizing  compressional  plasma 

[1,2] 

heating  by  an  imploding  liquid  metal  liner.  Use  of  a liquid 

(rather  than  solid)  liner  makes  possible  controlled  repetitive 
implosions,  in  which  most  of  the  energy  required  to  drive  the  liner 
is  recoverable.  This  reduces  the  required  value  of  Q,  the  ratio  of 
the  thermonuclear  yield  to  the  driving  energy  supplied,  and  also  the 
effective  dimensions  and  total  energy  needed  for  the  device,  since 
both  scale  as  Q^. 

One  obvious  problem  in  implementing  the  concept  of  plasma 
compression  by  cylindrical  liquid  liners  stems  from  the  possibility 
of  Rayleigh-Taylor  instabilities  at  radius  r = R (t),  the  interface 
between  liner  and  plasma.  When  the  liner  slows  and  then  rebounds, 
there  is  an  effective  gravity  g^^^  ~ R_  pointing  toward  the  center  of 
the  device,  i.e.,  from  the  liner  toward  the  plasma.  Interchange 
instabilities  localized  near  this  surface  can  be  expected  to  grow  with 
typical  rates 

V = (kgeff)^,  (1) 

where  k is  the  wavenumber. 

These  modes  can  be  stabilized,  however,  by  means  of  liner 
rotation.  Since  viscosity  is  negligible  (Reynolds  numbers  ^ 10^ 

are  typical  for  liquid  alkali  metals  with  scale  sizes  and  Implosion 
velocities  characteristic  of  current  experimental  configurations), 
angular  momentum  is  conserved.  Thus  each  annular  lamella  of  fluid  has 
an  azimuthal  velocity  v “ 1/R,  so  that  the  centrifugal  acceleration 
scales  as  R"^.  The  effective  gravity  now  has  the  form 

- K - v^/R.  (2) 

If  the  second  term  dominates  at  R ” R_,  there  and  the  Rayleigh- 

Taylor  modes  are  stabilized.  This  stabilization  mechanism  has  been 
Note;  ManuwriDt  submitted  February  22,  1978. 
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[5] 

demonstrated  experimentally  using  water  liners  and  shown  to  agree 
with  theoretical  predictions. 

Previous  theoretical  investigations  of  liner  Raylelgh-Taylor 
Instability  have  used  fairly  realistic  models  of  the  liner  motion 
assuming  the  fluid  to  be  incompressible.  Nevertheless,  non-zero 
compressibility  can  play  a significant  role  in  the  evolution  of  the 
perturbations.  If  the  growth  time  v"^  derived  for  an  incompressible 
fluid  is  shorter  than  the  time  required  for  sound  to  propagate  a dis- 
tance comparable  to  the  wavelength,  the  growth  rate  should  be  reduced. 
The  criterion  for  compressibility  to  be  important  is  thus  approxima- 
tely k ^ g/c^,  where  c is  the  sound  speed  in  the  liner  (typically  a 
few  km/sec  in  liquid  metals).  For  extremely  short  wavelengths,  on 
the  other  hand,  the  formula  [Eq.  (l)]  derived  for  an  incompressible 

basic  state  should  be  a good  approximation. 

[d] 

In  the  Linus-0  device^  and  in  any  larger  device  intended  to 
operate  in  a reactor  regime,  the  liner  undergoes  significant  compres- 
sion. This  is  particularly  true  near  the  inner  surface  at  turnaround. 

[7! 

Although  such  behavior  has  been  studied  in  detail  numerically,  no 

simple  analytic  model  of  the  basic  motion  is  available.  A linear 

r 4 1 

stability  analysis  analogous  to  that  of  Barcilon,  et  al.,  who 
described  the  evolution  of  the  perturbation  amplitudes  in  terms  of 
Integrals  over  the  dynamic  trajectories,  is  therefore  out  of  the 
question. 

In  this  paper  we  adopt  a simple  (and  consequently  somewhat 
unrealistic)  model  of  the  unperturbed  liner  motion  in  order  to  be  able 
to  get  exact  solutions  to  the  stability  problem.  This  model,  that  of 

isothermal  self-similar  contraction  and  expansion,  was  originally 

[bT 

derived  by  Korobeinikov  in  an  astrophysical  context.  In  common 

[pi 

with  a number  of  related  models  developed  by  Sedov  and  coworkers. 

It  yields  a time -dependent  nonlinear  solution  of  the  one-dimensional 
ideal  fluid  equations,  starting  with  the  assumption  that  in  Eulerian 
variables  the  velocity  u can  be  written  in  the  form 

u(r,t)  « rF(t).  (3) 
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The  density  p and  pressure  p are  then  determined  as  functions  of 

r^  and  quantities  derivable  from  F(t\  In  the  case  of  the  solution 

found  by  Korobeinikov,  p and  p are  connected  by  an  isothermal  law 

(speed  of  sound  c = c^ , independent  of  r and  t).  This  is  a good 

approximation  for  the  liquid  metals  in  Linus  devices,  where  the  maxl- 

mun  pressure  p ~ 10^^  dyne/cra^  ~ 0.1  Dc^  ^here  p is  the  density  at 
max 

zero  pressure),  and  consequently  the  corrections  to  the  lowest-nrder 
equation  of  state, 

p = ;P  - p), 

are  unimportant. 

The  basic  state  so  determined  describes  a nonrotating  liner.  Thus 
we  expect  to  find  instability  in  all  cases.  (If  sufficiently  rapid 
rotation  were  introduced  to  change  the  sign  of  g as  defined  in 
Eq.  (?),  stability  would  again  be  expected  to  occur.)  The  object  of 
the  present  work  is  to  see  in  a simple  model  how  compressibility  and 
the  time  dependence  of  the  unperturbed  motion  modify  the  characteris- 
tic dispersion  relation  Eq.  (P. 

A word  should  be  said  here  about  the  meaning  of  stability  as  we 
use  the  term.  In  analyzing  perturbations  about  time-varying  basic 
states,  one  in  general  finds  non-exponential  time  dependence.  We  have 
chosen  to  regard  a system  as  unstable  if  and  only  if  the  asymptotic 
development  of  the  perturbed  displacement  § is  faster  than  the 
asymptotic  expansion  of  the  system,  i.e.. 


lim  1 ?U)  I / RU)  ”,  (5) 

^ 00  “* 

where  both  5 and  R are  associated  with  a particular  fluid  element. 

This  definition  does  not  preclude  the  possibility  that  in  a "stable" 
system,!  ?(t)  j / R(t)  increases  for  a finite  time. 

The  paper  is  organized  as  follows.  In  Section  II  wo  obtain  the 
general  Sedov  similarity  solution  for  an  ideal  fluid  obeying  an 
adiabatic  law,  using  the  I.agrangian  derivation  previously  described  by 
Keller. We  then  specialize  to  the  isothermal  case  and  discuss  its 
applicability  to  liner  dynamics.  In  Section  III  we  develop  the  theory 


of  a linear  perturbation  about  an  arbitrary  time-varying  basic  state, 
again  in  Lagrangian  variables.  We  apply  this  theory  to  the 
Korobeinikov  isothermal  model  discussed  in  Section  II  and  Investigate 
the  stability  of  flute  modes  (axial  wavenumber  » 0,  azimuthal  mode 
number  m ^ O).  We  show  that  the  resulting  problem  can  be  completely 
solved  in  terms  of  confluent  hypergeometric  functions.  The  only 
modes  which  are  unstable  during  both  implosion  and  expansion  phases 
are  associated  with  perturbations  for  which  both  the  divergence  and 
curl  of  5 vanish.  (Most  of  the  details  of  the  calculations  described 
in  Sections  II  and  III  are  relegated  to  Appendices.)  We  conclude  the 
paper  with  a brief  discussion  of  our  results. 

II.  BASIC  STATE 

We  start  with  the  continuity,  force  and  state  equations  in 
Eulerlan  form: 


Ot 


+ 


V(p 


v) 


0; 


P 


+ Vp  = 05 


p = P(P). 


(6) 

(7) 

(8) 


We  assume  that  the  region  within  which  Eqs.  (6-8)  hold  is  bounded  by 
two  free  surfaces,  where  an  external  pressure  supplied  by  a massless 
fluid  balances  p.  Let  the  unperturbed  motion  of  a fluid  element,  the 
position  of  which  is  r^  at  t-0,  be  given  by  R = R(r^,t).  We  can 
then  transform  to  Lagrangian  variables,  i.e.,  write  the  displacement, 
density  and  radial  velocity  of  the  fluid  element  as  functions  of  r^ 
and  t.  The  similarity  hypothesis  becomes  in  Lagrangian  variables 
simply 

R - r^  f(t).  (9) 

We  assume  one-dimensional  motion  in  the  direction  of  symmetry,  the 
coordinate  of  which  is  denoted  by  R. 


We  find  (Appendix  A)  Chat  when  Eq,  (4)  Is  assumed  for  Eq.  (8)  and  we 
Introduce  dimensionless  variables,  density  being  measured  In  units  of 

A 

p,  velocity  In  units  of  c and  length  In  units  of  the  density  scale 
length,  Chen 


ff  - 1, 

(10) 

where  we  have  taken  f(0)  “ 1,  f(0)  “ 0,  and 

p(r^,t)  - (D/f)  exp  (-r^2/2)  - 1 + p(r^,t). 

(11) 

where  D Is  a constant.  Equations  ( 9~ll)  hold  within 
occupied  by  the  liner,  defined  according  to 

the  region 

(12) 

By  Eq.  ( 9 ) the  motion  of  the  two  surfaces  is  given  by 

R_  - R(r_,t)  - r_  f; 

(!’•) 

ilk) 

The  similarity  solution,  Eqs.  ( 9“11),  should  be  viewed  as  a 
convenient  model,  not  an  exact  description  of  the  motion  of  a real 

liner.  Nevertheless,  we  will  show  that  there  is  a parameter  regime 

In  which  the  model  is  entirely  appropriate,  namely,  in  the  immediate 

vicinity  of  the  turning  point  reached  by  a thin  liner  compressing  a 
plasma  or  magnetic  flux. 

Figure  1 displays  the  behavior  of  f as  a function  of  time. 

Note  that  f(-t)  ■ f(t),  so  that  the  trajectories  described  by  the 
similarity  model  are  symmetric  under  reflection  about  t » 0.  This 
contrast  with  the  previously  cited  results^^^  determined  by  solving 
Eqs.  (6-9)  self-consistently  with  an  adiabatic  gas  payload  and 
vanishing  external  pressure.  They  Indicated,  among  other  things,  that 
the  liner  expands  more  slowly  then  it  contracts,  due  to  tamping  by 
the  still-imploding  outer  layers.  The  latter  effect  is  least 
noticeable  for  small  peak  pressures,  1 e.,  when  p - p-1  « 1.  For  the 


5 


liner  trajectory  in  our  model  to  be  consistent  in  the  neighborhood 
of  the  turning  point  with  an  adiabatic  gas  law  for  the  plasma  being 
compressed,  viz., 

p - C(r./R_)2''  = Cf-^Y,  (15) 

we  must  have 

(o/f)  exp  (-r_^/2)  - 1.  { l6) 


Here  C Is  the  maximum  plasma  pressure  reached,  expressed  in  the  same 

A 

units  (Pc^)  as  the  similarity  solution  itself.  Using  the  early-time 
approximation  (A22)  for  f we  find  that 

D = (C  + 1)  exp  (r//2)  ( 17) 

and 

c = 1/(2y-1).  (18) 

For  a mixture  of  plasma  aid  magnetic  flux  it  has  been  shown^^^^  that 
the  pressure  is  well  approximated  by  a law  of  the  form  (l4)  with  an 
effective  value  of  V satisfying  5/5<Yg££  < 2,  whence  from  (l3) 

lA  < C < 5/7.  ( 19) 


The  strict  applicability  of  the  similarity  solution  as  a model 
for  the  liner  motion  is  evidently  limited  by  the  fact  that  for  suffi- 
ciently large  f the  pressure  becomes  negative.  This  occurs  first  at 
r * r^,  and  defines  a restriction  on  the  similarity  factor  to 

f * f , where 
max'^ 


fmax  * (C  -»■  1)  exp  [i(r_^  - r_^2)].  (20) 

For  thick  liners,  l.e.,  when  - r_^  >2^  (1+C);  Eq.  (20)  implies 
f <1,  which  is  Inconsistent  with  Eq.  (lO).  In  the  thin  liner 
limit,  when  - r ^ « 1,  the  model  can  only  be  strictly  applied 
for  times  such  that 

f s f « C + 1,  (21 ) 

max  ’ 


and  only  then  if  the  external  pres- >re  satisfies 
CgU)  - [(C  + l)/f]  exp  [i(r_^  - r^^)]  - 1 

(C  + l)/f  - 1.  (22) 


In  the  remainder  of  this  paper  we  will  relax  the 
restrictions  (19-22^.  Instead,  we  will  concentrate  on  finding  the 
complete  time  evolution  of  perturbations  about  the  basic  state,  in 
the  belief  that  their  behavior  will  display  the  general 
features  displayed  by  more  realistic  models.  Support  for  this  is 
provided  by  the  result,  established  below,  that  the  eigenfunctions 
of  the  most  unstable  modes  arc  sharply  peaked  at  r = r , i.e.,  far 
from  the  outer  layers  where  the  similarity  solution  is  unphysical. 

III.  PERTURBATIONS 


From  Appendix  B,  we  find  the  equation  in  Lagrangian  variables 
for  the  perturbed  displacement  §: 
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V + 


where  the  gradient  operators  on  the  right  hand  side  are  defined  with 
respect  to  the  instantaneous  unperturbed  displacement  R.  At  this 
point  we  assume  that  ^ is  independent  of  Z,  the  axial  component  of  R. 
As  a result  of  this  assumption,  there  is  no  fixed  spatial  scale  in 
the  problem  associated  with  a parallel  wavenumber.  Thus  the  gradient 
operator  V which  appears  on  the  right  hand  side  of  Eq.  (2'*i''  can  be 
written 

u 6^16 

^"^r6^-*-^eR  m 


+ ^0 


(2I4) 


We  can  then  multiply  Eq.  ^25^  by  f^,  leaving  an  expression  involving 
only  operations  with  respect  to  the  Lagrangian  position  (,  the  position 
of  C*0^  and  the  azimuthal  angle; 
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(25) 


f2  ^ - V (V  . ^ - r • £)  + i 

In  Eq.  (24)  and  the  remainder  of  Section  III,  we  suppress  the 
subscript  zero  In  r^.  The  initial  conditions 

i^(r)  = i(r,o),  {26) 

i^(r)  = i(r,0'  (27) 

are  arbitrary  except  that  they  must  satisfy  the  boundary  conditions 
which  hold  at  r = r^  , 

V . i = 0,  (28^ 

a consequence  of  assuming  that  the  external  pressure  is  unperturbed 
(Appendix  B).  The  problem  defined  by  Eqs.  i25-28)  can  now  be 
completely  solved. 

We  begin  by  rewriting  Eq.  (25^  in  the  form 

f ^ ^ + _r  • V ^ = W - r X ( 2'^ ) 

where 

a = V • ^ (70) 

and 

ui  = V X i.  (31) 

We  can  find  equations  for  a and  « by  operating  on  (29)  with  diver- 
gence and  curl: 

a = - a - r.  Vo+r.  Vxu)_;  (32) 

f^  uj  = ^.  (33) 

We  first  solve  Eq.  (35)  for  x.  Knowing  we  can  solve  Eq.  (32)  for 
<7;  using  the  results  found  for  ^ and  ct,  we  finally  solve  Eq.  i,20)  for 
5^.  The  details  of  the  analysl  s are  presented  in  Appendix  C.  The 
complete  solution  can  be  written  in  the  form 
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+ *0^  (£>t)  + e Y(r,t)  + (5I4) 


where 


a(r,0,t)  - i:  e 
m,n 

Im  u)  (0) 


im9 


n f*  n 

H +1 
r n 
ri 


"„(0) 


Im  iu„(0)  1 r im  ou  (0) 

-ir;r-  - " 

n J L 


M*  +1 
n 


, « (r) 

Im  n 

11  +1 
n 


P ( r,  6,  t ) = E e 
m,  n 

Im  (U  (0) 


Ar)  1 

f(t)  +i^(o)  g(t) 

li_  ^ drW  (r) 


lm9 


im  r 


I drw  ^ r 1 r 

" r — ^ \o 

“ M.  +2  ) I n 

r n < 
ra 


(0^ 


im  n r im  u)  (0)  T I 

■ir^]  Vo)-ir7r-  V'' 

n n ’ 

o ^ drW  (r)  , ^ 

-nr  !>  n .Ip 

r 1:  7 '> 


11  +1 
n 


drW^( r ) 


(“-(o)  f(t)  + g(t) 


(’i6) 


Y(r,0,t)  = E e 1 

m^  n f ' 


and 


«5(r,e,t)  - E e^"®  Ca_p  (t)  + B q (t)] 

m I " “ “ \r  / 


(^8) 


Here  W satisfies 
n 


“n"  * <i  - ')  - 1 - 7>"„  - ° 


{-^9) 


9 


and 


w„(r^)  ■ 0. 


where  [i  is 
n 

a (t) 

n 


V'' 


an  associated  eigenvalue; 

?rr  ^ 

= 5^]*  dee"^™®  drre"’^‘  W^{r)o^^,e,t); 

0 r_ 

2TT  ^ 

“ ^ dBe"^"*®  / drre  W (r^e  .u;^r,0,t); 

2n  <J  n -z  - ’ 

0 r- 


2TT  r_^ 

^ J dee”^"*®  J drre"’^  ,rW(r,e,t); 
0 r 


(i+O) 

(J*l) 

{k2) 


S (t)  = 5 [-i  i fi  t)]; 

n n 

T (t")  = f(t'l  4 [i  - i M'  , T ; ^ f{t')];  (41+) 

n n 

g(t)  = f(t)  / dse"®  (1+5)  • | 

0 ii 

I 

P^(t)  = 5 [i  (l^^m),  4;  f(t)l;  (1+6)  ^ 

q_^(t)  » f(t)  4 i;^f(t)],  (1+7) 


and  A^,  are  constants.  Here  4 Is  the  standard  solution  of  the  first 
kind  (Kunnner  function)  of  the  confluent  hypergeometrlc  equatlon^^^ 

All  of  the  quantities  defined  In  Eqs.  (kO-^7)  except  g depend  on  m, 
but  the  mode  number  label  has  been  suppressed  to  streamline  the 
notation. 

To  Investigate  stability,  we  consider  the  Implosion  and 
expansion  phases  separately.  First  we  note  that  f,  and  p^  are  even 
functions  of  t,  while  g,  and  q^  are  odd.  Solutions  of  pure  even 
(odd)  parity  are  symmetric  (antisymmetric)  about  t«0.  A linear 
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combination,  e.g.,  Q(Eq.(C5l)]  of  the  two  can  be  chosen  such  that 
the  solution  Is  a monotonlcally  Increasing  function  of  time 
asymptotically  for  both  t “ and  t "*  + ®. 

We  consider  first  the  outward  motion  starting  at  t=0.  From  the 
asymptotic  formula  (C5^)  we  see  that  the  only  time  dependence  In 
Eq-  (5^)  which  grows  faster  than  f(L'  Is  that  in  p (t^  and  q (t). 

[in  Appendix  D we  show  that  this  statement  is  true  for  all  t > 0, 

not  just  asymptotically.]  While  £ could  be  small  or  even  vanish  at 
t=0,  owing  to  cancellation  among  the  various  terms  in  the  solution, 
so  that  the  perturbations  could  grow  over  some  finite  time  interval, 
the  relative  amplitude  of  the  perturbations  vanishes  as  t ”*  * except 
for  modes  which  have  finite  A or  B in  Eq.  (58)*  These  are 
evidently  incompressible  irrotatlonal  modes  localized  at  the  inner 
surface. 

From  Eq.  (C54),  P_/f  and  q /f  both  grow  as  (£«  f')'"^^as  t 

For  t close  to  zero,  Eq.  (C25)  says  that  the  perturbations  grow 

exponentially  with  growth  rate  Y = (m-l'l^,  in  agreement  with  the 

[4] 

result  for  Rayleigh-Taylor  instability  in  an  incompressible  liner. 

Now  we  turn  to  the  question  of  the  .stabi lity  of  the  inward 
motion.  We  ask  under  what  circumstances  all  of  the  components  of 
^/f  can  have  smaller  values  at  some  t < 0 than  they  attain  at  turn- 
around. We  are  thus  led  to  seek  solutions  which  damp  with  increasing 
-t,  i.e.,  with  increasing  f.  Evidently  they  must  have 
and  y(0)  = “Y(0)  and  must  contain  the  functions  p and  q only  in  the 
combination  Q,  defined  in  Eq.  (C51)»  The  fastest  growing  modes  are 
those  with  all  constant  coefficients  set  to  zero  except  those  multi- 
plying S and  T . According  to  Eqs.  (C26-27^  and  (C54),  these  grow 
n n 

i.M.n  m+i+i  (|^'l 

as  f)®  ^ , or  using  Eq.  (C24),  as  {k/nf)  . Clearly 

the  amplification  Increases  with  both  n and  m.  These  modes  are  very 

different  from  the  familiar  Rayleigh-Taylor  instability.  They 

represent  overstable  internal  waves,  and  are  entirely  due  to  the 

presence  of  compressibility. 
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If  we  look  for  modes  which  are  unstable  both  before  and  after 
turnaround,  we  see  that  the  only  ones  are  those  with  Y * ^ = <7  = 0 
and  time  dependence  proportional  toQ(t).  Figure  2 aompares  the 
behavior  of  0*"  for  m » 1 and  10  with  that  of  f. 

IV.  CONCLUSION 

Using  Lagrangian  coordinates  throughout,  we  have  obtained 
general  formulations  of  the  one-dimensional  self-similar  motion  of 
an  ideal  isentropic  fluid  (already  in  the  literature^ ) and  of  the 
development  of  infinitesimal  perturbations  around  an  arbitrary  time- 
dependent  basic  state  (apparently  new).  The  two  are  combined  and 

applied  to  a treatment  of  the  Rayleigh-Taylor  instability  in  a self-  ll 

similar  fluid  motion.  When,  as  is  true  here  for  z-independent  per- 
turbations, the  gradient  operator  acting  on  the  perturbations  is 
homogeneous  in  r,  we  can  separate  the  space  and  time  dependence  of 

® great  simplification.  Thus  the  analysis  presented  here  of 
flute  perturbations  about  a cylindrical  self-similar  motion  is  act- 
ually simpler  than  the  corresponding  slab  problem,  which  involves  the 
perturbation  wavelength.  In  fact,  the  spherical  problem,  which  is 
treated  elsewhere, can  be  solved  for  all  modes  without  restriction 
on  the  angular  dependence  Y™(9,cp). 

For  liquid  metal  liners,  the  appropriate  relation  between 
pressure  and  density  is  Eq.  (J+),  the  isothermal  law.  Carrying  out 
the  analysis,  we  obtain  the  complete  solution  for  J^(r,t^,  Eq.  (^^). 

One  class  of  modes  grows  prior  to  turnaround  and  damps  afterward. 

These  are  internal  compressional  (sound^  waves,  which  are  pumped  up 
by  the  liner  contraction  and  relax  again  as  the  liner  expands.  The 
only  modes  which  grow  (relative  to  the  liner  radius''  both  before  and 
after  turnaround  are  those  which  are  Incompressible  and  irrotational. 

For  sufficiently  large  mode  numbers  m,  or  times  sufficiently  close  to 
turnaround,  the  growth  is  exponential  and  identical  to  that  found  for 
an  Incompressible  medium,  essentially  given  by  Eq.  (l).  Otherwise 
their  time  dependence,  given  by  Eq.  (C5l\  is  nonexponential.  Higher 
mode  numbers  m (shorter  wavelengths'*  give  rise  to  faster  growth,  with 
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the  amplification  relative  to  the  unperturbed  liner  radius  scaling 
asymptotically  for  large  t as  f)"*  ~ t)  . 

This  Is  similar  to  results  obtained  In  magnetohydrodynamlc 

r 121 

treatments  of  the  stability  of  static  equilibria.  There  it  is 

generally  found  that  the  fastest-growing  perturbations  are  those 
which  are  Incompressible,  since  compression  reduces  the  free  energy 
available  for  producing  magnetic  field  deformation  and  fluid  motion. 
It  is  also  Interesting  to  compare  our  result  with  the  well-known^  ^ 
solution  to  the  problem  of  the  Raylelgh-Taylor  instability  of  a 
stationary  isothermal  exponentially  stratified  medium  supported  from 
below  by  another  medium  of  much  lower  density,  where  it  is  found 
that  the  growth  rates  are  given  exactly  by  Eq.  (l).  This  conclusion, 
which  may  appear  surprising  in  light  of  our  arguments  (Section  l) 
that  compressibility  can  reduce  growth  rates  when  k ^ g/c^,  finds  an 
explanation  when  reached  by  the  techniques  of  Section  II. 

Equation  (28''  states  that  at  free  surfaces  the  medium  deforms 
without  compression.  Since  in  a static  exponentially  stratified 
medium  the  eigenfunctions  vary  exponentially  in  the  direction  normal 
to  the  boundary,  condition  (28)  holds  everywhere  if  it  holds  on  the 
boundaries.  In  effect,  it  is  propagated  inward  from  the  surfaces. 

In  the  problem  treated  here,  the  radial  dependence  of  the  perturba- 
tions is  not  exponential,  but  the  unstable  modes  are  essentially 
localized  at  the  inner  surface.  Hence  we  should  expect  to  find 
^ ” 0,  even  though  the  model  allows  compressible  motion. 
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APPENDIX  A.  SELF -SIMILAR  SOLUTIONS  OF  THE 
IDEAL  FLUID  EQUATIONS 


VIhen  we  introduce  Lagranglan  variables,  Eq.  (6)  can  be 

[12] 

Integrated  at  once  to  yield 


P(lo>0) 

“ 17  R(r  ,t)  1 ' 


(Al) 


while  Eq.  (7)  becomes 


P p dp  ^ 


1 do  < V> 


_ 1 ^ 


P dP  IV^I  • ^o  P ” - IV^l  • ''o  P dp 

(A2) 

where  dots  denote  derivatives  with  respect  to  t and  the  superscript  A 
denotes  the  adjoint.  Here  we  have  made  use  of 


( VL  .V  ; ^ ^ , 


V r • V R a I , 
~o  o—  = ' 


(A5) 


where  I^  is  the  unit  dyadic.  Hence 


7p  » [V  r^(R, t)]  • V^p  - 

R)^ 

R)‘"  • ^.P  • ^-P- 


0-'  o ‘iVl  ’ ° 


(aJ^) 


These  expressions  simplify  greatly  for  the  cases  of  planar,  cylindri- 
cal and  spherical  symmetry,  viz.. 


p dr  » pdR; 

O 0 

{A5^ 

2TTp  r dr  - S^pRdR; 

^ a6) 

O 0 0 

l^TTp  r ^dr  ■ ^TTpr^dRe 

0 0 o 

(A7) 

1 

i 


3^ 


15 


Equations  {A'j-'J)  can  be  summarized  In  the  single  relation 


p(r,t)  - p(r 


(a8) 


where  the  dimensionality  V = 1,  2,  3>  respectively. 

When  we  Introduce  the  similarity  assumption  Eq.  (9),  Eq.  (A8) 
reduces  to 

.-U 


P = P f 
o 


(A9) 


Moreover,  if  we  specialize  Eq.  (8)  to  an  adiabatic  law 


p = P [(p/Pr-1], 


(AlO'l 


where  p,  P and  Y are  constants,  then 


i2=V2/p\  ^2. 

dp  A I * I ^ A 

P \P/  P 


(All^ 


and  one  can  write  (A2)  in  the  form 


r f 
o 


dp  5 

dp  6r 


■Urn  (p/p)  = 


-»^(Y-l)  ^/P 

" f a 


or  if  Y ^ 1, 


. of  . . ^ ^ ^ 


£2 


Y-1 


(A12) 


P o 


where 


A A 


C2  « p/p 


(AlU) 


and  > 0 is  a separation  constant.  If  Y*1  (the  case  of  interest  in 
the  present  work),  Eq.  (AI3)  is  replaced  by 


ff  - c2/£2 


(A15) 
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Note  that  by  assuming  a positive  separation  constant,  ve  have 
effectively  chosen  case  III  (in  the  notation  of  ref.  9),  in  which  all 
fluid  particle  trajectories  pass  through  a minimum  displacement  at 
the  same  time  (when  f - O).  This  instant  we  denote  by  t-0. 

Equation  (Al^)  thus  yields 

0 £r>  p 

6(r^^/2)  ■ ■ (Al6) 

and 


1/f, 


(A17) 


where  we  have  adopted  units  such  that 

C-/-P-P-1.  (Ai8) 

The  solution  of  (Al6)  is 

Pq  " D exp  {~x^f2),  (A19) 

where  D is  a constant.  Equation  (A17)  can  be  integrated  once  to 
yield 


a 

f 


(2  Sn 


whence 


ds  e 


s2/2 


i 

r A 

J ds  e 


(A20) 


(A21) 


0 


0 


For  early  times  t « 1, 


f(t)  ••  1 + i t^. 


(A22) 


Integrating  (A21)  by  parts  repeatedly,  we  find 

f^/2 


t - ^ + — 


f-1  f 

O sees 

(2  ^ f)^ 


(A25) 


17 


A 


To  lowest  order,  then,  for  t 1 

___  r 1 T 

f(t)  ~ t t + 

(A2U) 

From  (A9)j 

p - p^/£  - exp  {-r  ^/2) 

o ^ 

{A25) 

and 

p ■ p-1  - exp  (-r  ®/2)  - 1. 

(A26) 

Finally,  the  radial  velocity  becomes 


APPENDIX  B.  DERIVATION  OF  THE  PERTURBED 
EQUATION  OF  MOTION 


Consider  the  equations  of  motion  In  Lagranglan  form. 


Eqs.  (Al-2).  Let  R(r^,t)  denote  the  unperturbed  motion  of  a fluid 


element  vhose  location  Is  r at  t«0.  Hence 

— o 


p^{^t)  R » " \ Pq, 


(Bl) 


vrhere  and  are  the  unperturbed  density  and  pressure,  and 


P(R,t) 


W R (r  ,l:)l 


{B2) 


Let  a {>erturbed  motion  be  defined  by 


(B5) 


Note  that  on  taking  the  gradient  of  (B5)  with  respect  to  r at  fixed 
t one  has 


^ = V R + y 

But  by  the  chain  rule 


(b4) 


V = (7R)  . = (I  - V i)  . 


■ ’r  - S'  • \ 


'r  - S'  ■ ’r 


(B5l 


Thus 


\\ l\  =•  l^R  R \ ll 


\l  + \i\  ^ I + \ ' I • 


(b6) 
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The  Integral  of  the  perturbed  continuity  equation  (the  condition  that 
mass  is  conserved  under  the  perturbation  B’>)  reads 

P(r,t)  - (B7) 

where  we  use  R to  identify  a given  Lagrangian  fluid  element. 

Together  (b6)  and  (.87^  yield 

P ‘R,  t) 

P„(R,t)  - D^{R,t^  V_  • 5 (b8) 


Hence  the  perturbed  density  is  given  by 


Now  the  linearized  equation  oi  motion  can  be  written 


V p 
R *^o 


(BIO'' 


i-  "o 'r  • s]  * ‘'r  5'  • • S'  V, 

["dp  dp  1 

R I dp  R ^ o ^ ^ R o do  I R ^ R *^o 

— o — — oJ  — — 


- ,v  • 5^  V p 
' R ^ R ^^o 


'r  \ ■ '»o  5']  - 'r  ts  • 'r  Xo’  * >'r  S'  • 'r  -c 


(V  . ?'!  V p 
' R R ^o 


I — • (p  “ 5 • p - • 5 p 

R [dP^  R '"o^J  ^ RR^^o  R ^rPc 


(Bll'l 


r 

If  we  divide  by  0^  and  rewrite  (Bll),  suppressing  the  subscripts 
R and  zero,  there  results 

S “ ^ [p  ^ ^ ^ S ^ 

- 1 • ^ ^p]  - ^ I • ^0  Vp  - 1 (V  . £)  7p 


Using  (A25)  and  (A26)  for  0 and  p,  we  find 

5-  [^.  .PS.]  .^S 

= V [v  . i . + S-  , U15) 

the  Lagrangian  form  of  the  perturbed  equation  of  motion  in  the  case 
of  an  isothermal  pressure  law. 

The  boundary  condition  at  a free  surface  associated  with 
Eq.  (B12)  is  found  from  requiring  that  the  pressure  on  the  boundary 
not  change  under  the  perturbations: 

p(R)  = ^(R)  + p^(R) 

■ PJS'  - 3?  'o'  • 5 1®“*' 

whence 

V . i « 0 (B15) 


APPENDIX  C.  SOLUTION  OF  THE  PERTURBATION  EQUATION 

The  Solution  of  Eqs.  (25-28''  is  facilitated  by  two  circumstances: 
the  space  and  time  dependence  separate  completely  (in  Lagrangian 
variables),  and  the  solutions  can  be  carried  out  in  stages  so  that  ^ 
is  found  "a  piece  at  a time."  The  first  piece  is  the  position 
vorticity  which  satisfies  Eq.  a homogeneous  equation.  One 

solution  is  Immediate,  since  from  Eq.  1,10), 

f^  f = f.  (Cl) 


The  other  is  then  found  to  be 

t t(t) 

..X  - I*  dt  ^ f . -s^/2 

0 0 

We  thus  have 

fiO)  = 1,  ho)  = 0; 

g(0^  - 0,  ho)  = 1, 


(C2) 

(C^) 

(CJ4) 


and  so  we  can  write  the  complete  solution  of  (55)  including  Initial 
conditions,  in  the  form 


i«(r,t''  = (iLgiz)  fit)  + ^ g(f)^ 


(C5) 


where  (“^(z)  = tu(r,0).  Evidently  the  spatial  dependence  of 
lu  is  "frozen  in"  at  t=0. 

We  now  Fourier  transform  (£  and  o in  angle  according  to 


2TT 


I™  (r,t)  = y a dS  e'^™®  i 


(c6) 


etc.  [in  what  follows  we  will  display  the  dependence  on  m in  the 
transformed  quantities  only  where  it  is  important.]  If  we  now  write 


i(r>t)  - e^a(r,t''  + Cq  e(r,t)  + e^  Y(r,t), 


.07) 
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subject  to 


W{r^)  - 0, 


(Cl6) 


determines  an  eigenvalue  problem  with  solutions  corresponding  to 


n ■ 1,  2,  . . ..  Writing  Eq.  (C15)  In  the  form 


1 d r -r2/2  dwl  ^ r , vF 
; 37  37j  * 


-r®/2 

e W - 0, 


(C17) 


we  see  that 


> 0,  (Cl8) 


dr  r e 


-t^/2 


W‘ 


so  that  the  eigenfunctions  are  oscillatory.  From  classical  Sturra- 


Llouvllle  theory,  we  know  that  the  eigenfunctions  W (r)  form  a 

n 


complete  orthonormallzable  basis  In  the  space  of  real  functions  on 
(r_,  r^),  with  the  orthogonality  condition  given  by 
r , 

_.-2  /o 

(C19) 


J dr  re"’^  W„(  r ) , ( r)  - 6 ^ , . 


n n 


nn 


The  W's  can  be  expressed  In  terms  of  confluent  hypergeometrlc 
functions,^ viz.. 


W^(r)  - r"*^  F Cl^  (1  ± m - p.^),  1 ± m;  i r^],  (C201 


where  F(a,  c;  x)  Is  a solution  of 


d^  F 


dx2  (c  - x)  g - a F 


0. 


(C21'' 


The  M.'s  must  in  general  be  determined  numerically.  The  most  efficient 
method  Is  to  utilize  (Cl8)  as  a variational  principle  and  solve  for 
by  the  Raylelgh-Rltz  procedure.  Given  p.,  a linear  combination  of 


1 
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1 


two  independent  solutions  of  (C2l)  can  be  found  which  satisfies  (C12). 
The  asymptotic  forms  of  the  standard  solutions  of  the  confluent 
hypergeometric  equation  (C2l)  in  the  limit  a "*  - * for  fixed  x 

,[15] 


are 


}(a,c;  x)  ~ n r(c)  (i  cx  - ax)  ^ exp  (ix) 
• cos  I [(2c  - i^a)  x]^  - i ttc  + i n I 


and 


i £ 


Y (a,c;  x)  ~ n F (-^  c-a  + ^)  x exp  (ix) 

• cos  I [ (2c  - Ua)  x]^  +(a-'ic  + 4)Trj 


(C22) 


(C25) 


Using  these  we  can  calculate  the  asymptotic  behavior  of  the  eigen- 
values M.^,  i.e.,  the  WKB  approximation  for  M'jj.  The  result  is 


n 


(C2U) 


n = 1,  2,  . . where  Ar  » r^  - r . All  but  the  first  few  eigen- 
values are  well  approximated  by  (C2U),  and  the  approximation 
Improves  as  Ar  “*  0. 

The  time  dependence  is  given  by 

f^S+^i„S-0.  (C25) 


Using  ^ f as  the  new  Independent  variable,  we  can  convert  (C25)  into 
a confluent  hypergeometric  equation  also.  The  result  is 

s(t)  - F(-  i f)  - f^  D (f),  (C26) 

n 

where  D is  a solution  of  the  parabolic  cylinder  (Weber-Hermite^ 
equation.  It  is  convenient  to  choose  two  independent  solutions  of 
(C26)  as 
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Sn{t)  - ^ ^ f),  (C27) 

T^(t)  - t(t)  ^ , i;  ^ f^,  (C28) 

where  t is  the  Kuimner  function^ associated  with  (C21).  From  this 
choice  It  follows  that 

S(0)  - 1,  S(0)  - 0;  (C29) 

T(0)  - 0,  T(0)  - 1.  (C50) 


We  next  seek  a particular  solution  of  the  Inhomogeneous  equation 
(Cll)  in  the  form 

a(r,t)  - Zo  (t)  W (r).  (C51'' 

n " " 

Using  (C5)  and  C15),  we  find 

f2  o + p,  a ■ Im  [u)  (O)  f(t)  + ti.  (O)  g(t)l.  (C52^ 

n n n n n 


The  desired  particular  solution  is  thus  simply 


a 

n 


Im 

p.  +1 

n 


[u)^(O)  f + g], 


(C55) 


and  the  complete  solution^  Including  Initial  conditions.  Is 

Im  (u  (O) 

l»  (0) 


c„(t)  - 


Im  (U  (O) 


s_(t)  + 


T (t)  + 77^  ClU.lO)  f(t) 


p +1  ■ n' 

n 


+ \(0)  g(t)l . 


(c^U) 
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^ Now  we  can  find  a particular  solution  of  (29)  in  the  form 
S,*ej.a+^  0,  where 


a(r,t)  - E r ” / 


drW  '(r) 


r^ 


c^n(O) 


im  U)  (0)  1 

-ir^l  V 


r.  im  i (0)  1 4 W (r) 

+ \o  (0) r T (t)  + -» 


tA  +1  r 
n 


(0)  - 


• j«»„(0)  f(t)  +<i^(0)  g(t)j  ; (C55) 

u ^ drW  ( r ) I r 

.<,,0- r I*.. » 0,- 

im  U)  (0)1  r im  u>  (o)"!  i 

V')-  [V°>-17^]  V'>| 

~ r J ~r f(tj  + iLjo)  g(t)l  + 

1 r / 

-S^  dr  rW^(r)  u,^(o)  f(t)  + <Ljo)  g(t)  (C56 

3 ' ' 


(C56) 


If  the  arbitrary  lower  limits  on  the  integrals  are  chosen  so  that 


<V*> 


7-5  ; (C37) 


+1  - r 2 


W ^ ( r ) * r W ( r ) , 

"3  an'  3''' 


(C58) 


(C39) 


u>  » 7 X 


(ck)) 


and 


V • I . 


(cui) 


Thus  the  general  solution  to  (29)  can  be  written 
S - g + 7 0, 


where 


72  0 - 0. 


th 


The  m component  of  0 satisfies 


whence 


1 & / 60”\  A® 

r Si  - 


■ *:  (tJ  - (r) " 


0 


(Ck2) 


(Cl+5) 


(CU) 


(cU5) 


The  time  dependence  Is  obtained  by  substituting  (C^5)  In  (C29), 
whereupon  the  right  hand  side  vanishes: 


f^  0 ® - ( 1 m)  0™  - 0 


(Ck6) 


This  equation  Is  of  the  same  form  as  (C2^),  so  we  can  Immediately 
write  down  the  solutions.  For  consistency  with  our  previous  choices, 
we  set 


P™(t)  - *Ci(l=Fm),  4;  ^ f], 

q^(t)  “ f 5 [1^  4m,  4;  ^ f]. 

These  satisfy 

P±(0)  - 1,  p±(0)  - 0; 

q±(0)  - 0,  q^(0)  - 1. 


(cl+T) 

(cJiB) 

(CU9) 

(C50) 
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For  m^k,  reduces  to  f multiplied  by  a polynomial  In  fin  f of  order 
k. 

Since  p_  and  q_  both  diverge  the  same  way  at  large  t.  It  la  use- 
ful to  Introduce  In  their  place  two  new  solutions  which  behave 
asymptotically  like  ▼,  the  standard  confluent  hypergeometrlc  function 
of  the  second  kind: 


P(t) 

Q(t) 


^ r(nri-l)  q (t) 


(C51) 


For  large  argument  we  have 

P(t)  ~ (fin  f)"  t - » ; (C52) 

Q(t)  ~ (fin  f)"  t -*  ® . (C55) 

All  the  other  asymptotic  results  we  require  follow  from  the  general 
formula 

S/  \ 1 z a-c 

(a,c;z)  ~ e z 


as  z "•  ® . 


(C54) 


i 
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APPENDIX  D.  A BOUND  ON  THE  FUNCTIONS  S . T 

n'  n 

Given  > 0,  let  U(t)  end  h(t)  satisfy 
u + n U - 0; 


h - h ■ 0, 

where  £ is  defined  in  Appendix  A,  and  let 

U(0)  - h(0)  - a > 0; 


U(0)  - h(0)  - b > 0. 


where  a + b - 1.  Then  the  following  inequalities  hold  for  0 < t < ": 


lu(t)l  S h(t)  S ^ f(t). 


Proof:  (1)  Assume  a-1,  b-0.  Then  h-f.  Denote  the  solution  of  (Dl) 
with  these  initial  conditions  by  U^.  Define 


gCu]  - ^ + u^  > 0. 


Since 


we  have 


-2^iU2f 


^ + U®  S liU2(0)  + U(0)®, 


or  for  U ■ U^, 


U^2  < f2. 


Hence  (D^)  holds  for  case  (1). 


Proof  supplied  by  Richard  Beals. 
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